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ABSTRACT 

We investigate the evolution of coUisionally merged stars with mass of ^ IOOMq which might be 
formed in dense star clusters. We assumed that massive stars with several tens Mq collide typically 
after ~ IMyr of the formation of the cluster and performed hydrodynamical simulations of several 
collision events. Our simulations show that after the collisions, merged stars have extended envelopes 
and their radii are larger than those in the thermal equilibrium states and that their interiors are 
He-rich because of the stellar evolution of the progenitor stars. We also found that if the mass-ratio 
of merging stars is far from unity, the interior of the merger product is not well mixed and the 
elemental abundance is not homogeneous. We then followed the evolution of these collision products 
by a one dimensional stellar evolution code. After an initial contraction on the Kelvin-Helmholtz 
(thermal adjustment) timescale 10'^"^ yr), the evolution of the merged stars traces that of single 
homogeneous stars with corresponding masses and abundances, while the initial contraction phase 
shows variations which depend on the mass ratio of the merged stars. We infer that, once runaway 
collisions have set in, subsequent collisions of the merged stars take place before mass loss by stellar 
winds becomes significant. Hence, stellar mass loss does not inhibit the formation of massive stars 
with mass of ~ lOOOM©. 

Subject headings: globular clusters : general — stars : early type — stars : evolution 



1. INTRODUCTION 

Recent infrared observations of the Galactic cen- 
ter and the centers of other nearby galaxies have 
revealed a population of compact and massive 
star clusters located close to the centers of galax- 
ies, such as the Arches and Quintuplet clusters 



(lOkuda et al.l [l990l : 
iFieer et all 120021 ) . IRS 
and 



iFiger. McLean fc Morris 



13E (Maillard et al 



and 



Lu et al.l [20051) 

1 12001 . 



MGG-11 in M82 
The estimated 
the range of 10'* Af© 
radii are between 0.1 
densities in excess of 



IRS 16SW 

(jMcCradv. Gilbert, fc Graham 
masses of these clusters are in 
to 1O^M0 while their half-mass 
and 1 pc, giving rise to central 
IO^Mq/pc^. 

Dynamical simulations have also shown that if 
star clusters are born with sufRciently high central 
density, massive stars with > 20Mq will sink to 
the cluster center within a few Myr, i.e., before the 
end of the stable nucl ear burning phase, through 
dyna mical fri ction ([P ortegies Zwart fc McMillan I 
2002t iPortegies Z wart ct al] 



2004 



Freitag. Rasio fc Baumgardt 12006 
Freitag. Gtirkan fc Rasio, .2006! ) . Their stellar radii 
are large enough that there is a high chance for collisions 
between them after the stars arrived in the center. 

Indeed, iV-body and Monte Carlo simulations have 
shown that collisions between high-mass stars in young 
star clusters can lead to the formation of a supermassive 
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star with mass of several hundreds to several thousands 
M„ ("Portegies Zwart et al.''2 004l : lBaumgardt et al.ll2006t 
[Freitag. Gtirkan fc Rasio 2006[). Such supermassive stars 
and the intermediate-mass black holes (IMBHs) which 
might form out of them could be the ultra-luminous X- 
ray sources re cently discovered by Chandra and HST 
observations (Hopman. Portcgics Zwart. fc Alex ander I 
[2004; Baumgardt et al. 2006; Patruno ct al. 2006). It is 
also argued that runaway massive stars that are ejected 
from d ense clusters become a candidate of gamma-ray 
bursts ([Hammer et al.ll20Q6f ). 

However, whether supermassive stars can really form 
through runaway collisions and whether IMBHs forms 
at the end of their lifetime is still not clear. So far, 
most simulations have neglected hydrodynamical pro- 
cesses during the collisions and the effects of stellar 
evolution. Stars formed from the merging of other 
stars might start their lives with significant abundance 
gradients because of incomplete mixing. In addition, 
since merging events happ en only after a star cluster 
has gone into core-collapse (IPortegies Zwart et al.l 120041 : 
iFreitag. Gtirkan fc Rasiol |2006[ ). runaway stars initially 
have higher He abundance as a result of the nuclear burn- 
ing of the parent stars. Evolution of stars formed through 
merging, thus, is likely to be different from the evolution 
of homogeneous stars with "normal" abundance. 

It has been argued that the formation of IMBHs 
from metal-rich (~ solar abundance) massive 
stars is u nlikely , because stron g ste llar winds 
([Kudritzki I l2002t iNugis fc Lamers I |2000D consider- 
ably reduce the masses before th e blac k holes form 
jBelkus. v an Bever fc VanbeverenI 120071 : lYungelsonl 
2006). However, because the lifetime of merged stars 
with higher He content is shorter than that of normal 
stars, the total mass lost during the lifetime might be 
smaller than these estimates. 
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The present paper is a first attempt towards a realis- 
tic treatment of the stellar evolution of runaway stars. 
In the present paper we will follow the collision of two 
stars by means of SPH calculations and then follow the 
evolution of the merger product with a stellar evolution 
code. The parameters of the colliding stars are taken 
from the results of iV-body simula tions of runaway merg- 
ing of stars in young star clus ters (jPortegies Zwart et al.l 
[200I iBaWgardt et al.|[2006l ). 

2. METHOD 

Our procedure consists of three steps : (1) stellar evo- 
lution of single stars, (2) simulations of stellar collisions 
and (3) stellar evolution of collision product s. This pro- 
cedure is es sentially same as that used by iSills et aL~l 
()1997l I2OOID . In their works, they concentrated on the 
formation of blue straggler stars due to the collision be- 
tween low mass stars (< IMq). However, in the present 
paper we are interested in the merging process and sub- 
sequent stellar evolution in the core of a very dense star 
cluster where only massive stars (> lOM©) are involved 
in collisi ons because of mass se gregation. 

jPortegies Zwart et al.l (|2004l ) showed that collisions in 
very dense clusters typically start at t ~ IMyr after their 
formation. Following this result, we consider various col- 
lisions (Table [T|) of massive stars after t = 1 Myr from 
the formation. We determine the interior structure of 
the merged stars by calculating the stellar evolution of 
single stars with solar abundances from zero-age main 
sequence (ZAMS) phase (Step (1) of the procedure). 

The stellar evolution is handled by one-dimensional 
(ID) spherical symmetrical stellar evolution code. Our 
code, which adopts a usual Henyey metho d, is based 
on th e program originally developed by iPaczvnski I 
()197Cl( ). As adopted in general stellar evolution calcula- 
tion, we neglect hydrodynamical evolution and only treat 
evolution on a Kelvin-Hclmholtz (thermal adjustment) 
timescale, 

where G is the gravitational constant, and M^, (Mq), 
Ri, (Rq), and L^, (Lq) are stellar (solar) mass, ra- 
dius, and lumin osity, respectively. Our code adopts the 
OPAL o pacity (llglesias &: Rogersi [T9 96) and equation 
of state ( Rogers. Swenson. &: Iglesias 1996) tables, and 
the nuclear reaction rates tabulated in Bahcall & Ulrich 
(|T988( ). We focus on the early phase of stellar evolution 
before the central H is exhausted. For this purpose it is 
reasonable to switch off nuclear burning of He and heav- 
ier elements, and to consider only H-burning because the 
central temperature is still not high (< 5 x lO^K in our 
simulations). 

For Step (2), we construct spherical symmetric three 
dimensional (3D) stellar structures from the ID results. 
Then, we perform 3D smoothed-particle-hydrodynamical 
(SPH) simulations of stellar collisions with the parame- 
ters summarized in Table [U We use a mod i fied y ersion 
of the SPH code by .Nakasato fc Nomoto I (|2003[ ). An 
important modification is that we include a treatment 



of equation of state (EOS) that takes into account ra- 
diation pressure in addition to a fully ionized ideal gas 
(7 = 5/3). We adopt a Balsara type artificial viscosity 
(jBalsara III995D with viscosity parame ters a = P = 5/6 
as suggested by Lombard! et al. (2003) . 

For the equal mass cases (EQl and EQ2 in Table[T]), we 
use N = 10, 000 particles for each star (mass resolution 
of 8.85 X 10~^ Mq). For the unequal mass cases (UEl 
and UE2), we use N = 20,000 particles for Star 1 and 
N — 6305 for Star 2, respectively (mass resolution of 
4.42 X 10^"^ Mq). As will be discussed in Appendix this 
relatively small N is sufficient for our purpose and our 
results of SPH simulations do not strongly depend on N. 

We assume that initially each star is separated by im- 
pact parameter, 2 x (i?i + R2), in all the cases, where 
i?i and i?2 are the stellar radii before the collision. We 
put the first star at the origin and the other at the x-axis 
with a specific tangential velocity as shown in Table [TJ 
These velocities are typical for c ollisions that occur in 
centr al regions of dense clusters (jPortegies Zwart et al.l 
l2004f) . Larger initial velocities lead to less eccentric or- 
bits and larger pericenter distances between the two stars 
and therefore longer merging times. During SPH simu- 
lations, we always check whether merging occurs or not 
by assigning each particle into a membership either 1st 
star, 2nd star or unboun d using enthalpy of each particle 
(jRasio fc Shapiro I fl99ll ). We stop a SPH run at least 
^imorging = 2 x 10^ scc after the mering. Note that the 
choice of ^imcrging IS arbitrary and does not affect the re- 
sults as long as it is sufficiently longer than the dynamical 
time scale. This criterion ensures that the particles are 
in dynamically stable states. The last snapshot is used 
to create the structure of a merged star that is evolved 
in Step (3). 

In Step (3), we follow the evolution of the collision 
products by the same ID stellar evolution code used in 
Step (1). To do this, we need to construct ID structure 
from the results of 3D SPH simulations in Step (2). In 
this paper, we neglect the effect of rotation, and sim- 
ply average radius, r, density, p, temperature, T, and 
elemental abundances of the 3D results in mass radius 
coordinate, m to give ID spherical symmetric structure. 
Needless to say, this treatment is a very simplified one. 
The SPH simulations show that the end-products of the 
merging events are not spherical due to the rotation (see 
Figure [6]). For example, the ratio between the final rota- 
tional velocity and the circular velocity at a given radius 
ranges between 0.2 (inner region) and 0.8 (outer region). 
Indeed, for more detailed studies, we should adopt a more 
elaborated way that takes into account rotation when 
mapping a 3D d istribution into a ID spherical symmet- 
ric profile (e.g.. ISills et al. ||2001D . However, since the 
rotational energy of the merged stars is less than 15 % 
of the gravitational energy, the rotation is not expected 
to affect the evolution of the stars much. Thus, we think 
that we can give rough but reasonable estimates for the 
evolution of the collision products by our simple prescrip- 
tion. 

On the other hand, the SPH simulations cannot treat 
low density envelopes on account of limited mass reso- 
lution. Then, we extrapolate density and temperature 
structure to the outer region to match the inner struc- 
ture obtained by the SPH simulations. In order to do 
this, we adopt the Eddington approximation to derive 
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EQl 


EQ2 


UEl 


UE2 


Star 1{Mq) 


88.5 


88.5 


88.5 


88.5 


Star 2{Mq) 


88.5 


88.5 


27.9 


27.9 


e 


0.444 


0.100 


0.669 


0.125 


yinit(kms-l) 


550 


700 


400 


650 


collision time (day) 


2.1 


44 


1.4 


13 


Final Mass (Mq) 


165.5 


156.4 


106.1 


98.0 


-dM{MQ) (collision) 


11.5 


20.6 


10.3 


18.4 



TABLE 1 

Summary of the four runs. The third and fourth lines 

SHOW the orbital eccentricity e AND INITIAL TANGENTIAL 
VELOCITY OF THE STARS. THE LAST LINE GIVES THE MASS LOST 
DURING THE COLLISIONS {—dM). 



the relation between temperature and optical depth and 
assume hydrostatic equilibrium to set the density struc- 
ture in the outer envelopes. 

Also, in the SPH simulation we do not consider nuclear 
energy release, because the total nuclear energy produc- 
tion integrated with collision duration is much smaller 
than the gravitational energy. In the later evolution of 
merged stars (Step 3), however, nuclear burning is as 
important as the energy release by gravitational contrac- 
tion. Then, we have to carefully determine luminosity, I, 
from an energy equation. 



dl_ 
dm 



du\ p 



dp 
dt 



du \ dT 



dT J dt' 



(2) 



where Cn is net energy gain by nuclear burning minus 
neutrino loss and u is internal energy (e.g. section 9 
of Kippenhahn & Weigert 1990). Note that the terms 
involving time derivative denote the energy release (ab- 
sorb) by gravitation contraction (expansion) . As a mat- 
ter of calculation technique, we have to set an appro- 
priate time-step, At, for these terms in order to pro- 
ceed the stellar evolution calculation stably. We adopt 
At « (1/100)tkh(= 10 - lOOyr) (but much larger than 
free-fall time scale) just after the merger events to pre- 
cisely follow Kelvin-Helmholtz contraction. As merged 
stars settle down to thermally stable state, we adopt 
larger At(< tkh)- 

Because collision products are not in thermally equilib- 
rium states just after the mergers, we need to prepare an 
appropriate initial guess for the time evolution of phys- 
ical variables (p, T, and r as functions of m) by the 
Henyey method; without this treatment we fail to fol- 
low the evolution of collision products, even though an 
appropriate At is set. For an initial guess of the correc- 
tion, we use a mixture of the structure of SPH simula- 
tion and thermally relaxed structure (e.q. ZAMS) that 
is determined separately. Then, we derive the correct 
time-evolved structure by relaxation. 

3. RESULTS 

3.1. SPH simulation oj stellar collisions 

In all cases reported in the present work, the overall 
evolution is qualitatively similar: After a certain time, 
both stars are elongated by tidal interaction to form an 
extended merging product as shown in Figures [1] and [21 
which depict the last phases of the merging processes 
for cases EQl and UEl respectively. However, the mass 
ratio and the orbital eccentricity of the merging stars 
influence the details of the merging processes. 









Fig. 1. — The last phase of the merging process in the SPH sim- 
ulations for EQl case. In these snapshots, particles are projected 
onto the xy plane (orbital plane). The size of each panel is 100 
Rq. Blue and red points represent the particles originating from 
Star 1 and 2, respectively. 



3.1.1. Mass Ratio 

The mass ratio affects the material mixing in the in- 
terior of the merger products. Figure [1] shows that the 
material is well mixed in the equal mass case, and the 
elemental abundance is almost homogeneous as will be 
shown later. On the other hand, in the unequal mass case 
(Figure [2]) the less massive star (Star 2), which has the 
higher central density, sinks to the center without suffi- 
cient mixing, and the more massive star (Star 1) forms 
an extended envelope. This is consi stent wit h recent re- 
sults for lower mass stars bv lDale fc Davies I (i2006'). As 
a result, the elemental abundance is also inhomogeneous 
in the merger product. This difference of the material 
mixing affects the later evolution of the merged stars, 
which will be discussed in §3.21 

Let us examine the material mixing during the mergers 
in more detail. Figures [3] and [3] respectively compare the 
structure of colhsion products (soHd lines) of EQl and 
UEl with the corresponding parent stars just before the 
mergers (dashed and dotted lines). The ID structure of 
the merged stars are reconstructed from the 3D results 
at the ends of the step (2) as explained in ^ From top 
to bottom, density, p, pressure, p, an entropy variable, 
pI p"' ^ and H and He abundances are plotted against mass 
radius, m/MQ^ where 7 is a ratio of specific heats that is 
determined from the OPAL equation of state table. An 
entropy variable, p/p'*', is useful to study the material 
mixing; under the adiabatic condition convective stabil- 
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Fig. 2. — Same as Figure [T] but for UEl case. 

ity in uniform media with gravity in — z direction reads 
d{p/ p^) / dz > (Schwartzschild criterion), and in con- 
vection zone d{p/ p^)/dz = 0. Thus, a fluid element with 
smaller p/ tends to sink into a central region during 
merger. 

The third panel of Figure [3] shows that in the parent 
star (Star 1 — Star 2) of EQl the large convective core 
(constant p/ p^) extends to BOM©. In the He-rich enve- 
lope, p/p'^ is larger by only r:! 0.2 dex except the region 
that is very near the surface. Hence, the matter is well 
mixed during the merger because of the low entropy bar- 
rier. As a result of the mixing, the gradients of the H and 
He abundances are very small (the bottom panel Figure 
^. 

This result is in contrast to the collision o f low mass 
mam sequence (MS) stars. iLombardi eFan (|1995D per- 
formed the collision of two O.SM© MS stars. They found 
that the mixing between the dense cores and the outer 
envelopes was inefficient. This is because the interior 
of a low mass star is occupied by a radiative core and 
the fraction of a surface convection zone is tiny in mass. 
Thus, p/p'^ monotonically increases in the interior; e.g. 
in O.SMq, p/ p"^ increases by nearl y an order of magnitud e 
from the center to the envelope (jLombardi et al.ll2002f ). 
As a result, the material mixing is inhibited during the 
merging by the large entropy barrier, which gives a clear 
contrast to our result for massive stars. 

Next, let us move on to the collision of the unequal 
mass stars (UEl; Figure H]). The third panel shows that 
p/p'^ of the lower mass parent (Star 2; dotted line) is 
smaller than that of the higher mass partner (Star 1; 
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50 100 150 

m/MQ 

Fig. 3. — Stellar structure of the parent stars (dashed; Star 
l=Star 2 in this case) and the collision product of EQl at the 
end of step 2 (solid). The ID structure of the merged stars are 
derived from the 3D SPH results by averaging in m (see [(Jjl. From 
top to bottom, density, pressure, entropy variable, and H and He 
abundances are plotted as functions of mass radius, m/MQ. 



dashed line) because the lower mass star has higher cen- 
tral density (the top panel). Owing to the small en- 
tropy the lower mass star (Star 2) sinks to the core of 
the merged stars without mixing, which was seen in Fig- 
ured] As a result, the elemental abundances in the core 
reflect those of Star 2, while the abundances in the outer 
region reflect those of Star 1 (the bottom panel); the 
merged star consists of the H-rich core and the He-rich 
envelope, which is opposed to usual evolved single stars. 

3.1.2. Orbttal Eccentricity 
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Fig. 4. — The same as Figure[3]but for UEl. The structure of Star 
1 and Star 2 are shown in dashed and dotted lines, respectively. 

The orbital eccentricity controls the time for the stars 
to merge. Figure [5] shows the evolution of the distance 
between the centers of both stars. In the run with ini- 
tially more eccentric orbit (initially smaller pericenter 
distance), the stars merge more quickly than those in 
more circular orbit, mainly because the kinetic energy of 
the system is smaller; EQl takes ~ 2.1 days to merge 
whereas EQ2 takes ~ 44 days. EQl loses only 11.5 Mq 
during the merging owing to the shorter collision dura- 
tion, while EQ2 loses a larger mass of 20.6 Mq. 

Despite the different merging time scales, the struc- 
tures of the merged stars, EQl and EQ2, are not so dif- 
ferent after they settle down to dynamically stable states 
(Figure ini). Therefore, the later evolution of the merged 
stars is not expected to be different between EQl and 



EQl 
EQ2 
UEl 

/ \ A s UE2 

I' !!;" >■ 
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Fig. 5. — Evolution of the distance between the centers of the 
merging stars. In EQl and UEl runs, the two stars merge after a 
few orbital revolutions, whereas the EQ2 and UE2 runs need much 
longer time to merge. 




Fig. 6. — Density snapshots of the merger product at the end 
of the SPH simulations. Particles are projected onto the xz plane 
where z is the rotation axis. The size of the panels is 62 Rq. Left 
: EQl and Right EQ2 

EQ2. We found similar tendencies in the unequal mass 
cases, UEl & UE2; in UEl (larger e), the merger product 
settles down to a dynamically stable state faster and the 
mass lost during the merger is smaller, while the later 
evolution is essentially the same. 

3.1.3. Rotation, Difference etc. 

As noted previously, the structures of the merged stars, 
EQl and EQ2 (also UEl and UE2), are rather similar but 
angular momentum (AM) distributions of the merged 
stars are slightly different. Figure [7] shows AM distri- 
bution as a function of mass radius for EQl and EQ2. 
Clearly, outer AM distribution of EQ2 is larger than that 
of EQl because of the difference in the merging time 
scale. Namely, during longer orbital revolutions in EQ2, 
more orbital AM is transfered to outer material than in 
EQl. Similar trends are observed in UEl and UE2 such 
that UE2 has slightly larger AM in that outer region 
than UEl. Although we neglect the effect of rotation 
in Step (3) (evolution of merging stars) in the present 
work, we expect this small difference in outer AM dis- 
tribution gives little influence on subsequent evolution of 
the merged stars. 

3.2. Evolution of Collision Products 

As explained earlier, we follow the evolution of the 
merged stars by the ID stellar evolution code. Since 
the mass ratio of the merging stars affects the structure 
of the merger products, and since the orbital eccentricity 
does not, we mainly study the evolution of EQl and UEl, 
and only briefly mention the results of EQ2 and UE2 for 
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Fig. 7. — Angular momentum distribution as a function of mass 
radius for EQl (solid line) and EQ2 (dashed line). 
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Fig. 8. — Evolution of stellar structure of EQl. The left panels 
show temperature (top) and density (bottom) as functions of radius 
in units of Rq . The right panels show entropy variable, %>l . (top), 
and H and He abundances (bottom) as functions of mass radius, m, 
in units of Mq. The dashed, dot-dashed (only for the left panels), 
dotted, and solid lines are the results at t = 0, 54 yrs, 6460 yrs, 
and 1 Myr after the merging, respectively. Note that the radius of 
the star is minimum at t = 6460 yrs. 



comparison. First, we study in detail the stellar evolu- 
tion without taking into account the effect of mass loss 
by stellar winds. Later in i i3.2.31 we present the results 
with mass loss for comparison. 

3.2.1. Equal Mass Collision 

Figure [8] shows the evolution of the stellar structure 
after the merging of the equal mass stars (EQl). On the 
left, temperature (top) and density (bottom) are plotted 
against t/Rq, and on the right entropy variable p/p'^ 
(top) and elemental abundances are plotted as functions 
of to/Mq. Just after the merger the envelope is ex- 
tended to i?* = 296Oi?0, exhibiting a core- halo struc- 
ture owing to the opacity peak around T ~ 2 x lO^K 
(jlshii. Ueno. fc Kato I Il999f) as shown in the left pan- 
els of Figure [S] The bottom right panel shows that the 
interior is well mixed during the merger; the elemental 
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Fig. 9. — Luminosity normalized by Lq as a function of m/MQ 
for run EQl. The solid and dotted lines give the total luminosity 
and the luminosity due to the nuclear burning at t = 0. At this 
stage most of the energy released comes from the gravitational 
contraction of the star. The dashed line shows the total luminosity 
at t = 6460 yrs, where the energy is coming from nuclear burning 
at this time. 



abundance is almost homogeneous, 0.52 < X < 0.6 and 
0.46 <Y < 0.38, even at t = 0, where X, Y are H and He 
abundances. This is mainly because the entropy gradient 
is small in the parents stars (^i j3.1.ip . 

Figure [9] presents luminosity, 1/Lq, as a function of 
m/MQ. The solid and dashed lines are the total lumi- 
nosity derived from Equation ([2]) at i = and 6460 yrs. 
The dotted line is the luminosity due to only nuclear re- 
action at t = 0. When we calculate this, we integrate the 
nuclear energy term (en) on the right hand side of Equa- 
tion Because the central density and temperature 
are lower than those in the thermal equilibrium state at 
i = 0, the contribution from the nuclear burning to the 
total luminosity is very small. Instead, most of the en- 
ergy comes from the gravitational contraction and it is 
transported outward by convection in the large convec- 
tive core up to mass radius, m ~ IIOM0; the star is in 
a state similar to that of a pre-main sequence star. The 
total Z at t = (solid line) decreases outward between 
mass radii 20 < m/MQ < 90. This is because the liber- 
ated energy is not converted to radiation but to internal 
energy, namely the increase of the temperature in this 
region. 

After the merging, the star contracts towards the ther- 
mal equilibrium state on a Kelvin-Helmholtz timescale 
(eq. [1]). The temperature and density increase (the left 
panels of Figure [8]), and the convective core, which is the 
region with constant p/p'^, also grows until it occupies 
« 80% of the total mass (the top-right panel of Figure 
[8]). As a result, the elemental abundances become homo- 
geneous inside m < I3OM0 (the bottom-right panel of 
Figure[5]). At i = 6460 yr after the merging, the star con- 
tracts to a minimum radius, i?^ = 37.6i?0, after which it 
expands gradually. At this time the star is in the stable 
H-burning phase, and the luminosity is all from the nu- 
clear reaction (Figure [9]). The later evolution traces the 
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Fig. 10. — Evolution of the radii of the merged stars, EQl 
and UEl (solid lines), in comparison with the chemically homoge- 
neous single star with mass, 165.6Af0, and abundance, {X,Y) = 
(0.6,0.38), and the star with IOG.IMq and (X, Y) = (0.61,0.37) 
(dashed lines). 




Fig. 11. — Evolution of the merged stars in a HR diagram. In 
the left panel, the solid lines show the results of the four cases. 
Stars are the locations just after the merging. In the right panel 
the evolution of the ch emi cally homogeneous single stars that are 
the same as in Figure 1101 (dotted lines) in comparison with the 
results of EQl and UEl (solid lines). The triangles are the initial 
locations (ZAMS) of the homogeneous stars. 

evolution of the chemically homogeneous single star with 
the corresponding mass and initial abundance, which we 
discuss below. 

Figure [TO] summarizes the evolution of the radii of the 
merged stars (EQl as well as UEl; UEl will be dis- 
cussed in §3.2.2p . in comparison with those of the sin- 
gle stars with the same masses and similar abundances 
= (0.6,0.38) for EQl and (0.61,0.37) for UEl). 
The figure clearly illustrates the monotonical contraction 
in t < 6460 yrs, which is followed by the gradual expan- 
sion during MS phase of stable H burning. The evolution 
during the MS phase is similar to that of the single star. 
The increase of the radii at i > IMyr indicates the end 
of MS phase due to the exhaustion of H. Our calculation 
shows that the MS lifetime of EQl is « 1.3 Myr, which 
is shorter than the corresponding lifetime (« 2Myr) of 
a solar abundance star with the same mass due to the 
smaller initial H abundance. 

Figure [11] shows the evolution of the merged stars in a 




R/Rg m/Mg 

Fig. 12. — Evolution of stellar structure of UEl. While cor- 
responding to Figure [8] for EQl, this figure focuses on the initial 
thermal adjustment phase. The solid, dashed, dotted, and dot- 
dashed lines are the results at t = 0, 34 yr, 149 yr, and 7490 yr, 
respectively, where in the top-right panel the dotted line is omit- 
ted to avoid confusion. Followed by the initial contraction before 
t < 34 yr, the star slightly expands between 34 yr < t < 149 yr 
and again contracts between 149yr < t < 7490 yr to the minimum 
radius. 

HR (Hertzsprung-Russel) diagram; the left panel shows 
the evolution of the four merged stars and the right panel 
compares the evolution of EQl and UEl with the corre- 
sponding single stars (for the unequal mass cases, see 
§3.2.2p . Through the initial contraction, the effective 
temperature, T^g, of EQl increases and L^, decreases. 
The position of the end of the Kelvin-Helmholtz con- 
traction (the turning point) in the HR diagram is near 
the ZAMS of the corresponding single stars (triangles). 
The later evolution resembles that of the single stars as 
we stated previously. 

The evolution of EQ2 (smaller e than EQl) is essen- 
tially similar to that of EQl (the left panel of Figure 
[TT|) . Because the mass is slightly smaller, the luminosity 
becomes lower (Figure [Tl]) . Note that Toff of EQ2 just 
after the stable nuclear burning sets in (the turning point 
in the HR diagram) is slightly higher than that of EQl, 
although the mass of EQ2 is smaller. This is because 
more massive stars (EQl) have more extreme core- halo 
structure to show lower T^s (Ishii et al. 1999). Chemi- 
cally homogeneous ZAMS stars of which masses exceed 
a certain limit also have this inverse trend. According 
to Ishii et al. (1999), Tcff of solar metallicity stars de- 
creases on increasing mass in the range of stellar mass 
> IOOM0, while less massive stars show the usual trend 
of the positive correlation between Tcs and stellar mass. 

3.2.2. Unequal Mass Collision 

The evolution of the merger products of the unequal 
mass stars is different from that of the equal mass cases 
during the initial contraction phase, while the later evo- 
lution follows chemically homogeneous single stars with 
the corresponding masses and abundances in both cases. 
Figures [T2l and [T3l present the results of UEl, which cor- 
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Fig. 13. — Luminosity normalized by Lq on m/MQ of UEl. 
The dotted line is the luminosity due to the nuclear burning at 
t = 0. The solid and dashed lines are the total luminosity at t = 
and 7490 yr, respectively, whereas the luminosity at t = 7490 yr is 
from the nuclear energy. 



respond to Figures [5] and [S] for EQl. As we have shown in 
the most important difference is that the interior is 
not well- mixed (Figures [2] and |4]). Just after the collision, 
the lower mass parent star (Star 2) sinks to the center 
without sufficient mixing because it has a higher density 
and lower entropy. This star is more H-rich since the nu- 
clear burning proceeds more slowly than in the massive 
partner (the bottom panel of Figure |4]). Therefore, the 
merged star consists of a H-rich core and a He-rich outer 
region as shown in the bottom right panel of Figure 1121 
Reflecting the higher density in Star 2, the density of 
the core of the collision product becomes slightly larger 
than the thermal equilibrium value, while the lower den- 
sity envelope extends to the outer region (the bottom- left 
panel of Figure [T2|) . 

Due to the high density as well as the moderate tem- 
perature, the nuclear burning takes place rather rapidly 
even just after the merging event in the unequal mass 
case (Figure [T3l) . The energy release rate by the nuclear 
reaction exceeds the total luminosity because the nuclear 
energy is also used to increase the temperature (internal 
energy; the top- left panel of Figure[T2|) and to expand the 
core (work on gas, i.e., the decrease of the core density; 
see the bottom- left panel of Figure [T2|) . Accordingly, the 
envelope also expands from t = 34 yr to 149 yr (Figures 
[TOl and fT2|l . In fact, a.t t = 149 yr, the radius becomes, 
i?i 412i?0, which is larger than i?^, = 337i?0 just after 
the merger (t — 0). Reflecting the initial expansion, the 
evolutionary path in the HR diagram (Figure [TT|l is also 
more complicated compared to the equal mass case. 

During this phase, the size of the convective core is 
small, m < 35Mq, because the entropy (the top-right 
panel of Figure [T2|) stays low in the core. Therefore, 
the H-rich core is still preserved without mixing with the 
outer region (the bottom right panel of Figure [T^ . In- 



cidently, the entire region outside m — 17 Mq becomes 
convectively stable with respect to the Schwartzschild 
criterion, d{p/ p"^) /dr > (the top-right panel of fT2)) . 
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Fig. 14. — Evolution of the stellar radii of EQl and UEl that 
take into account mass in the stellar evolution (solid lines), in com- 
parison with the results without mass loss (dotted lines). 



Between 17 Mq < m < 35Mq, the gradient of mean 
molecular weight, /i, leads to mixing because heavier He 
is more abundant in the upper layer; this region is un- 
stable only by the Ledoux criterion. Note that this is 
opposed to usual situations, in which heavier elements 
are more abundant in a lower region and /x gradient con- 
tributes to stabilization. 

The initial expansion between 34 < i < 149 yr is 
followed by the usual contraction to the equilibrium 
state through thermal adjustment. The chemical abun- 
dance becomes homogeneous from inward as the convec- 
tive core grows to to ~ 80Mq. The minimum radius, 
i?i — 2O.li?0, occurs at i = 7490 yr, roughly correspond- 
ing to tkh- The later evolution traces the evolution of 
the single homogeneous star, which is the same as in the 
equal mass case. The duration of the MS (« 1.6Myr) is 
again shorter than the corresponding lifetime (« 2.5Myr) 
for a solar abundance star. 

The evolution of UE2 is similar to the evolution of UEl 
(FigurefTTjl: The merger product initially consists of a H- 
rich core and a He-rich envelope. Although this structure 
is maintained at first due to the small convective core, 
the interior becomes homogeneous after t > 5000 yr as 
the convective core grows. The later evolution resembles 
the evolution of the corresponding single homogeneous 
star. 

3.2.3. Mass Loss by Stellar Wtnds 

So far we have not considere d the effect of mass loss by 
radiation-driven stellar winds (jCastor. Abbott, fc Klein I 
|1975[ ). However, it is supposed to affect the evolution 
of the merged massive stars. We study the evolution 
of EQl and UEl by explicitly taking into account mass 
loss in the stellar evolution calculations. Here we have 
adopted the mass loss rate M of solar metallicity gas 
from lKudritzki I ()2002[ ). which tabulates M as a function 
of Li,, Tcff, and metallicity^. 

^ Although in IKudritzki I l l200l) the depend ence of M on 
He abundance, Y, is not explicitly presented, observation of 
Wolf-Rayet stars show s that M has a dependence on oc yl-^-^ 
IINugis fc Lamers 120001) . Our merger products are He-rich Y 0.4 
in the envelopes, compared to the Sun (Y = 0.28), hence, M could 
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Fig. 15. — Evolution of EQl and UEl with (solid lines) and 
without (dotted lines) mass loss in a HR diagram. 

In Figures [T3] & [15] we show the evolution of EQl 
and UEl with mass loss, in comparison with the results 
without mass loss. Our collision products have typically 
M ~ (1 - 3) X 10-5 Mo/yr, and, ~ 10 - 30Mq is lost 
during the MS phase of ~l-2Myr. The differences of the 
stellar radii between the cases with mass loss and the 
cases without mass loss are not large (Figure [T^ ; they 
are less than 10% except at the very end of the MS phase 
(t > IMyr in EQl and t > 1.5Myr in UEl). Figure [H] 
illustrates that the luminosity becomes slightly smaller 
by ^ O.ldex at later epochs because of the mass loss. 

Once runaway collisions start in a dense cluster, the 
timesc ale of subsequent collisions is much shorter than 
IMyr (jPortegies Zwart et"alll2004f) . Therefore, further 
stellar collisions of the merged stars would take place 
before the mass loss becomes important in the stellar 
evolution. Metal-poor stars give even smaller mass loss 
rates than solar abundance stars. Thus, we can conclude 
that stellar mass loss does not stop the increase in mass 
due to runaway collisions, provided that the metallicity 
is comparable to or smaller than the solar value. 

4. CONCLUSION AND PROSPECT 

Bearing in mind formation of supermassive stars and 
IMBHs in dense star clusters, we have studied the hy- 
drodynamical processes during collisions of massive stars 
and the evolution of the merger products. After the colli- 
sions, the merged stars settle down to dynamically stable 
states on typical timescales of days to weeks, well before 
they would undergo further collisions. During the merger 
events, the stars typically lose ~ 10% of the total mass. 
The merger products are He-rich because of the nuclear 
burning of their parent stars. The interior of the merger 
product of equal-mass progenitors is well-mixed during 
this dynamical phase because of the low entropy barrier. 
On the other hand, during the merging of unequal-mass 
stars, the less massive star sinks into the core, and the 
more massive partner is elongated by tidal interaction to 

be be larger by a factor of 1.5-2 than that by IKudritzki I II2002I ). 
However, even if M becomes larger by this extent, we suppose that 
the effect of the mass loss is not still crucial during the MS phase. 



finally form the envelope. Since the nuclear burning took 
place slower in the less massive progenitor, the merged 
star consists of an H-rich core and an He-rich envelope. 

After the merger phase, the merged stars evolve 
to thermal equilibrium states on Kelvin-Hclmholtz 
timescales, lO^"'* yr. The evolution of the collision prod- 
uct of equal mass stars is very similar to a pre-main se- 
quence star; the star monotonically contracts and the 
luminosity is mainly supplied from the release of grav- 
itational energy. On the other hand, the evolution of 
the merger product of unequal mass stars is rather com- 
plicated due to the poorly mixed interior; the nuclear 
burning is already switched on owing to the sufficiently 
dense core, and as a result, the star slightly expands at 
first, which is followed by the usual contraction. 

Just after the merging the radius is larger than the 
equilibrium value by a factor of 10-100, while it goes 
down to a few times the equilibrium value in less than 
1000 years (See figure fTO]). N-body simulations that as- 
sume the mass-radius relation for MS stars show that 
collisions typically take place every ~ 3 x 10'*yr in a 
very dense region (Baumgardt et al. 2006). The col- 
lision probability will be enhanced when we take into 
account such realistic stellar radii. Then, a small frac- 
tion of merged stars might experience further collisions 
during the initial contraction phase if they are in a very 
dense region. However, most of the collisions are off-axis, 
and we expect that such off- axis collisions simply blow 
away a tiny fraction of the outer envelope, rather than 
resulting in a merger, because its density is quite low 
(Figures [5] and [TH). 

After the thermal adjustment phase, the merged stars 
enter a stable nuclear burning phase and their evo- 
lution is well approximated by those of single homo- 
geneous stars with corresponding masses and abun- 
dances. An important point here is that the lifetimes of 
merger products are shorter than solar abundance stars 
with the same masses because they are alre a dy H e-rich 
from the beginni n g. [Portegies Zwart et all (j2004f ) and 
iBaumgardt et all ( 2006h assumed 3 Mvr as the lifetime 
of runaway stars for the stellar dynamics simulations. 
Our result have shown that massive stars typically col- 
lide after IMyr from cluster formation and that the life- 
time of the merged stars is ~2Myr; these results confirm 
that the assumption of 3Myr is quite reasonable. This 
is robust even if merged stars experience further colli- 
sions because the collision products become more He (or 
heavier element)-rich and their lifetime is short. 

Our simulations show that neither mass loss during 
stellar collisions nor mass loss by the stellar winds pre- 
vents the growth in mass of the collision products. We 
can therefore anticipate that the scenario of the for- 
mation of supermassi v e stars by successive collisi ons 
(|Ebisuzaki et all 120011 : iPortegies Zwart et all l2004[ ) is 
likely to occur in realistic situations. Because of the 
nuclear burning, the merged stars become more H- 
poor. Finally, the lifetimes of massive He-rich descen- 
dants are much shorter than those of solar abundance 
stars with corresponding masses. We speculate that 
the very end-products of runaway collisions would form 
IMBHs quickly before suffering substantial mass loss; 
the key is that the material which is finally taken in 
a supermassive star spends most of the time in less 
massive stars, which are not influenced by mass loss 
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Fig. 16. — Structures of the merged stars EQl (left) and UEl (right) obtained with SPH simulations using different N. From top to 
bottom, density, internal energy and H abundance are plotted as functions of mass radius, m/MQ. Red and green lines correspond to 
results for = 20, 000 and 60, 000 in EQl and = 26, 305 and 78, 915 in UEl, respectively 



SO much. Therefore, supermassive stars are possibly 
formed by successive colhsions, ahhough it seems difS- 
cult through the evolution of very mas sive single stars 
(jBelkus. van Bever fc Vanbeverenll2007f l. However, our 
present work does not quantitatively treat this final pro- 
cess. For such purpose we need to the study evolution of 
very massive (~ IGOGM©) and chemically evolved (He- 
rich with abundance gradient) stars. We plan to carry 
out such simulations in the future. 
The gravity calculation of the SPH simulations 
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ence, and Technology of Japan. HB acknowledges sup- 
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APPENDIX 

A. DEPENDENCE ON NUMERICAL RESOLUTION 

In the present work, we choose rather small numbers o f parti cles {N 10,000 — 20,000) in Step (2). In previous 
extensive simulations of colliding stars bv lFreitag fc Beiizl (j2005[ ). they have repeated the same runs with different N's 
to see how the numerical resolution affects the hydrodynamical simulations of colliding stars. They have concluded 
that N ^ 10, 000 — 30, 000 per star is adequate to construct comprehensive tables for simulations of dense stellar 
systems. Our main interest in the present work, which is different from theirs, is the fate of collision products of 
massive stars as a result of the stellar evolution. 

The results of the resolution tests for EQl and UEl runs are shown in Figure [TBI All the quantities are averaged 
over every 55 particles in mass coordinate. In both cases, the runs with three times more N shows almost identical 
results to the corresponding original runs. Moreover, we practically average the quantities of SPH simulations of Step 
(2) to give ID spherical symmetric structure used in Step (3). Then, fine details in the results are smoothed out and 
the tiny differences seen in Figure [16] have little impact on calculations in Step (3) . We conclude that the relatively 
small N that we are using in the present work is sufficient for our purpose. 
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